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Abstract 

We  are  collecting  and  analyzing  large  and  diverse  datasets  of  seismic  waveforms,  absolute 
and  relative  body  wave  traveltimes,  and  surface  wave  phase  velocities  and  amplitudes.  Our 
primary  objective  is  to  map  the  variations  of  elastic  mantle  properties  beneath  Eurasia  over 
horizontal  lengthscales  of  approximately  1000-1500  km  and  vertical  lengthscales  of  50-100 
km.  In  comparison  with  recent  Harvard  global  models,  such  as  S12WM13,  we  have  increased 
the  horizontal  resolution  by  expanding  the  global  structure  up  to  spherical  harmonic  degree 
I  =  20,  and  we  have  adopted  a  radial  parameterization  which  gives  better  vertical  resolution 
in  the  upper  mantle. 

Using  earlier  and  new  datasets,  we  have  derived  a  new  whole  mantle  model,  S20U7L5. 
The  new  model  shows  many  of  the  same  features  seen  in  previous  tomographic  images,  but 
with  significantly  sharper  definition,  in  particular  in  the  upper  mantle.  Lithospheric  roots 
beneath  the  stable  continental  interiors  are  the  most  prominent  feature  in  the  uppermost 
mantle.  The  model  is  partially  validated  by  comparisons  with  the  earlier  regional  models 
SNA  and  TNA.  While  there  is  striking  similarity  in  the  velocity  structure  beneath  a  great 
variety  of  continental  interiors.  West  Africa  appears  to  overlie  one  of  the  fastest  locations  in 
the  upper  mantle,  with  faster  than  average  velocities  extending  to  500  km  depth. 

Additional  tasks  to  be  completed  before  the  new  model  is  finalized  is  the  incorporation  of 
more  accurate  corrections  for  crustal  thickness  and  structure,  as  well  as  the  inclusion  of  an 
extensive  dataset  of  recently  measured  teleseismic  traveltimes. 
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Objective 

We  are  adapting  the  seismic  tomographic  techniques  developed  in  global  scale  studies  to 
the  regional  scale  problem  of  mapping  the  elastic  and  anelastic  material  properties  beneath 
Eurasia.  The  main  objective  of  this  research  is  to  obtain  a  tomographic  seismic  velocity 
model  for  the  structure  beneath  Eurasia  with  a  horizontal  resolution  corresponding  to  at 
least  I  =  20,  consistent  with  the  widest  possible  range  of  seismological  observations.  Related 
objectives  are  to  evaluate  the  resulting  model  and  its  utility  in  event  location  algorithms. 

Research  Accomplished 

In  earlier  reports  (1993,  1994)  we  have  presented  results  based  on  earlier  whole  mantle 
and  preliminary  upper  mantle  models,  while  also  addressing  the  question  of  teleseismic  event 
locations  using  heterogeneous  mantle  models.  In  this  report  we  describe  in  greater  detail  the 
derivation  of  the  new  global  mantle  model  S20U7L5.  We  discuss  and  compare  some  shallow 
mantel  structures  seen  in  this  model  with  previously  derived  regional  models,  and  point  out 
some  remaining  problems/issues  which  we  plan  to  address  before  producing  the  final  version 
of  the  model. 

A  degree  20  model  of  the  Mantle 

Earlier  tomographic  models  of  the  whole  mantle,  based  on  a  variety  of  datasets,  have 
provided  descriptions  of  large,  planetary-scale  variations  of  P  and  S  velocities  in  the  mantle. 
With  the  goal  of  obtaining  better  resolution  of  structures,  particularly  for  the  upper  man¬ 
tle,  we  have  collected  and  incorporated  high-quality  surface  wave  dispersion  data  with  the 
previously  developed,  and  now  expanded,  seismological  datasets  of  travel  times  and  body 
and  mantle  wave  waveforms  in  the  model  inversion.  The  new  model  is  parameterized  to 
preferentially  resolve  smaller-scale  structures  in  the  upper  mantle,  while  at  the  same  time 
accounting  for  the  dominantly  large-scale  structures  in  the  lower  mantle. 

Parameterizati  o  n 

The  parameterization  and  inversion  approach  are  similar  to  those  used  for  earlier  Harvard 
mantle  models,  such  as  S12WM13  (Su  et  ai,  1993),  with  the  following  notable  differences: 

•  We  use  a  split  parameterization  for  the  mantle,  and  parameterize  the  radial  variations 
in  the  upper  mantle  in  terms  of  Chebyshev  polynomials  of  degree  0-7.  For  the  lower 
mantle  we  use  Chebyshev  polynomials  of  degree  0-5.  This  split  parameterization,  with 
more  degrees  of  freedom  in  the  upper  mantle,  reflects  our  goal  of  resolving  primarily  the 
shallow  mantle  structure,  and  the  greater  sensitivity  of  our  data  sets  to  heterogeneity 
in  the  upper  mantle  (see  below). 

•  The  spherical  harmonic  expansion  is  extended  from  degree  12  to  degree  20.  If  we  as¬ 
sume  that  we  can  resolve  structures  of  dimensions  that  are  approximately  half  the  min¬ 
imum  wavelength  of  the  highest  harmonic,  this  implies  a  resolving  length  of  ~1000  km. 
The  new  parameterization  increases  the  number  of  parameters  to  be  determined  for 
each  radial  degree  from  169  to  441,  and  the  total  number  of  unknown  coefficients  is 
14  X  441  =  6174.  While  we  do  not  expect  to  be  able  to  resolve  degree  20  structure 
in  the  deeper  mantle,  due  to  our  data  coverage,  we  do  believe  that  this  resolution  is 
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achievable  in  the  top  several  hundred  kilometers  of  the  upper  mantle.  We  damp  the 
spatial  derivatives  of  the  model  in  the  inversions,  and  a  higher  damping  is  applied  to 
the  lower  mantle. 

•  The  data  are  corrected  for  crustal  structure  using  a  new  compilation  of  crustal  thick¬ 
nesses  by  Mooney  (1994).  This  is  in  contrast  to  previous  models  which  were  only 
corrected  for  average  continental  and  oceanic  crustal  structure.  The  crustal  correction 
is  important  for  explaining  phase  velocities  of  shorter  period  surface  waves. 


Data  Sources 

Several  types  of  data  with  sensitivity  to  velocity  variations  in  different  parts  of  the  mantle 
are  included.  So  far,  the  following  data  have  been  incorporated  in  the  derivation  of  the  new 
model,  which  we  call  S20U7L5  (spherical  harmonic  degree  20,  degree  7  in  the  upper  mantle 
and  degree  5  in  the  lower  mantle): 

•  Absolute  and  differential  travel  times  measured  from  waveform  data  (S,  SS,  ScS, 
SS—S).  These  are  the  same  data  sets  that  were  used  in  the  development  of  S12WM13. 

•  Long-period  body  and  mantle  waveforms  from  451  earthquakes  recorded  on  the  Global 
Seismic  Network.  This  dataset  is  similar  in  coverage  to  the  data  that  were  used  for 
S12WM13,  but  we  use  a  more  broadband  spectral  weighting  in  the  waveform  fitting. 

•  Dispersion  measurements  of  Rayleigh  waves  at  140,  100,  and  75  seconds  (described 
below). 

•  Dispersion  measurements  of  Love  waves  at  140,  100,  75,  60,  and  50  seconds  (described 
below) . 

The  Love  wave  and  Rayleigh  wave  dispersion  datasets  have  not  previously  been  used  in 
inversions  for  mantle  structure,  and  they  are  described  in  some  detail  in  the  next  section. 

Dispersion  of  Love  and  Rayleigh  waves 

The  method  developed  by  Ekstrom  et  al.  (1993)  matches  an  observed  surface  waves  seis¬ 
mogram  o{t)  with  a  model  seismograms  m{t)  through  iterative  adjustment  of  parameters 
which  describe  the  dispersion  and  amplitude  spectrum  of  the  model  waveform.  A  fundamen¬ 
tal  mode  model  waveform  m{t)  is  calculated  using  the  JWKB  expressions  for  surface  waves, 
which  reduce  to 

m{u)  =  A(a;)  exp[z$(a;)].  (1) 

The  propagation  part  of  the  phase  is  expressed  as 


ujX 

c°{u)  -I-  6c{oj) 


(2) 


where  ^c(a;)  is  the  average  phase  velocity  perturbation  with  respect  to  the  PREM  value 
c°.  The  amplitude  function  A{uj)  and  the  pheise  velocity  perturbation  Sc{u)  are  parame¬ 
terized  as  smooth  B-spline  functions  of  the  period  band  of  interest,  35-150  seconds.  Since 
the  adjustment  of  m{t)  to  resemble  o{t)  cannot  be  determined  by  direct  linear  inversion,  an 
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iterative,  non-linear  algorithm  was  developed.  The  algorithm  is  based  on  two  concepts:  Iter¬ 
ative  Frequency  Band  Expansion,  and  Iterative  Phase  Isolation  and  Minimization  of  Residual 
Dispersion,  which  are  integrated  in  the  following  algorithm: 

1.  Spline  parameters  for  surface  wave  dispersion  and  amplitudes  are  calculated  for  an 
initial  Earth  model,  for  example  S12WM13,  and  a  trial  model  seismogram  is  calculated. 

2.  Observed  and  model  seismograms  are  filtered  between  75-125  sec,  the  initial  frequency 
window. 

3.  Dispersion  and  amplitudes  are  adjusted  to  maximize  the  fit  in  the  current  frequency 
window. 

•  Observed  and  model  seismograms  are  cross  correlated,  and  the  model  seismogram 
is  auto  correlated. 

•  An  objective  function  is  formed  which  reflects  the  difference  between  the  cross- 
and  autocorrelations. 

•  Spline  coefficients  for  5c{u})  and  A{u))  are  iteratively  adjusted  to  minimize  the 
objective  function. 

4.  The  frequency  window  is  expanded  to  higher  and  lower  frequencies. 

5.  Steps  3  and  4  are  repeated  until  the  whole  frequency  band  of  interest  is  included. 

Ekstrom  et  al.  (1993)  automated  this  algorithm  and  applied  it  systematically  to  surface 
waves  recorded  on  GSN,  CDSN,  GEOSCOPE,  and  MedNet  stations.  Moderate  and  large 
earthquakes  over  a  4-year-period  have  now  been  used,  and  more  than  80,000  paths  were 
analyzed,  out  of  which  approximately  10,000  Rayleigh  and  8,000  Love  waves  were  of  high 
quality.  Ekstrom  et  al.  have  reported  on  the  derivation  from  these  data  of  phase  velocity 
maps,  but  in  the  derivation  of  S20U7L5,  the  raw  phase  anomalies  were  used  and  treated  in 
the  same  way  as  travel  times. 

Inversion  Results 

The  new  model  S20U7L5  shows  many  of  the  same  features  seen  in  previous  tomographic 
images,  but  with  significantly  sharper  definition.  Figure  1  shows  a  comparison  between  the 
S  velocity  perturbations  at  150  km  depth  beneath  the  Mediterranean  and  North  Africa  as 
imaged  in  S12WM13  (top)  and  in  the  new  model  S20U7L5  (bottom).  Many  features,  such 
as  the  very  fast  West  African  craton  and  the  slow  velocities  coincident  with  the  Cameroon 
line  are  much  clearer  in  the  new  model.  The  slow  velocities  associated  with  the  Red  Sea 
rift  are  also  much  better  defined  in  S20U7L5.  While  the  new  model  shows  much  greater 
lateral  contrasts  in  the  upper  mantle,  the  low  angular  degree  components  of  the  new  model 
in  general  agree  very  well  with  previous  mantle  models  such  as  S12WM13. 

On  a  global  scale,  the  fast  velocities  beneath  stable  continental  interiors  is  the  most  promi¬ 
nent  anomaly  in  the  top  300  km  of  the  mantle.  The  very  high  velocities  beneath  the  stable 
interior  of  West  Africa  continue  even  deeper,  as  can  be  seen  in  a  velocity  profile  (Figure  2). 
The  new  model  also  shows  that  the  slow  velocity  anomalies  associated  with  ridges  continue 
to  great  depth,  in  particular  in  the  Indian  Ocean.  The  Red  Sea  rift  shows  slow  velocities 


620 


extending  down  to  400  km  (Figure  2).  Below  250  km,  however,  there  are  other  slow  features, 
not  directly  associated  with  surface  tectonics,  which  predominate,  in  particular  in  the  center 
of  the  Pacific  Ocean. 

Even  though  the  tomographic  models  result  from  a  formal  inversion,  damping  and  incom¬ 
plete  knowledge  of  data  and  model  uncertainties  make  it  difficult  to  assess  the  true  resolution 
of  our  models.  For  example,  a  question  which  might  be  asked  is  whether  the  deep  extent  of 
fast  velocities  beneath  West  Africa  is  well  resolved  and  imaged.  If  it  is,  the  velocity  struc¬ 
ture  beneath  this  shield  is  different  from  that  seen  beneath  the  Canadian  shield  (Grand  and 
Helmberger,  1984),  where  faster  than  average  velocities  appear  to  terminate  at  200-250  km 
depth. 

One  test  which  can  be  made  is  to  see  how  Canadian  shield  structure  is  imaged  in  S20U7L5. 
Figure  3  shows  a  comparison  between  the  vertical  velocity  profile  SNA  (Grand  and  Helm¬ 
berger,  1984),  representative  for  the  Canadian  shield,  and  S20U7L5  evaluated  at  a  point 
just  South  of  Hudson  Bay.  The  agreement  between  the  two  is  striking,  both  in  terms  of 
the  overall  amplitude  of  the  perturbation,  and  the  depth  at  which  heterogeneity  becomes 
small.  Note  that  the  velocity  jump  at  220  km  is  in  our  starting  model  (PREM),  and  that 
the  smooth  perturbations  appear  to  counteract  it.  Also,  there  is  agreement  between  SNA 
and  S20U7L5  on  the  vertical  gradient  below  400  km,  which  disagrees  with  PREM. 

A  similar  comparison  can  be  made  between  the  model  TNA  (Grand  and  Helmberger, 
1984),  and  S20U7L5  (Figure  3).  The  agreement  is  not  as  striking  as  for  the  shield  structure, 
which  may  in  part  be  due  to  the  fact  that  TNA  was  derived  from  data  which  sampled  a 
variety  of  tectonic  environments,  including  the  East  Pacific  Rise,  in  addition  to  western  US. 

While  the  lithospheric  root  beneath  West  Africa  extends  deeper  than  that  beneath  the 
Canadian  shield,  the  velocity  structure  beneath  several  other  stable  continental  interiors  are 
very  similar  to  each  other,  as  well  as  to  SNA.  Figures  4  shows  a  composite  of  six  different 
velocity  profiles  which  show  similarities,  not  only  in  the  upper  few  hundred  kilometers,  but 
also  in  the  deepest  part  of  the  upper  mantle. 

Remaining  Issues 

The  model  S20U7L5  predicts  all  the  datasets  that  went  into  its  derivation  well,  indicating 
that  the  observations  are  compatible  and  complementary.  Three  concerns  remain  to  be 
addressed  before  we  finalize  the  model: 

•  Are  the  corrections  for  shallow  structure  (crustal  thickness,  etc.)  adequate?  While  lin¬ 
earized  corrections  are  adequate  for  longer  period  data  (see  Woodhouse  and  Dziewon- 
ski,  1984),  we  will  investigate  the  effects  of  using  more  exact  crustal  corrections  for  the 
intermediate  period  surface  waves. 

•  Is  the  lower  upper  mantle  (250-600  km)  well  resolved?  We  have  recently  developed  a 
new  dataset  of  very  long-period  surface  waves  which  we  will  include  in  the  inversion. 
This  should  provide  better  constraints  on  the  structure  in  this  depth  range. 

•  Are  we  damping  the  structure  too  much  in  the  lower  mantle?  We  will  incorporate 
the  additional  datasets  of  all  major  teleseismic  arrivals  which  we  have  been  measuring 
from  waveform  data  over  the  last  two  years.  Several  of  these  phases  have  sensitivity 
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primarily  in  the  lower  mantle,  and  incorporation  of  them  will  improve  the  resolution 
there. 


Conclusions  and  Recommendations 

Intermediate  period  surface  waves  provide  excellent  global  coverage  and  sampling  of  shal¬ 
low  upper  mantle  structure.  The  inclusion  of  dispersion  data  for  these  waves  in  the  derivation 
of  the  new  model  S20U7L5  has  resulted  in  a  higher  resolution  and  better  constrained  model 
of  the  shear  wave  velocity  in  the  upper  mantle.  Several  interesting  features  are  seen  in 
the  new  model,  including  striking  similarities  of  structure  beneath  continental  cratons.  In 
West  Africa,  however,  fast  velocities  appear  to  extend  to  500  km  depth,  while  beneath  the 
Canadian  shield  the  fast  structure  terminates  at  around  300  km.  Ridges  also  show  a  great 
variability  in  their  expression  at  depth,  with  the  Indian  and  Antarctic  ridges  having  a  distinct 
expression  down  to  at  least  200  km  depth. 

Lateral  variations  in  S  velocity  in  the  upper  200  km  of  the  mantle  (as  imaged  in  S20U7L5) 
are  large,  -5%  -  -1-7%.  It  is  likely  that  in  order  to  resolve  heterogeneities  on  shorter  length- 
scales  than  those  seen  in  S20U7L5,  new  tools  have  to  be  developed  which  take  better  account 
of  wave  propagation  eflfects  in  heterogeneous  media.  In  addition,  improved  corrections  for 
crustal  thickness  and  crustal  structure  will  make  it  possible  to  better  resolve  mantle  veloc¬ 
ities  just  below  the  Moho.  Incorporation  of  additional  seismic  data  for  regional  distances 
(e.g.  Pn,  Sn)  should  also  help  to  improve  the  resolution  in  the  shallowest  mantle. 

By  combining  a  large  number  of  diverse  seismological  observations,  we  have  made  progress 
towards  obtaining  higher  resolution  tomographic  images  of  the  regional  scale  elastic  struc¬ 
ture  globally  and  beneath  Eurasia.  Large  and  regional  scale  heterogeneity  in  the  mantle  has 
been  shown  to  bias  seismically  derived  event  locations,  and  it  appears  likely  that  further  im¬ 
provements  in  locations  can  be  achieved  by  correcting  for  better  known  elastic  heterogeneity 
in  the  Earth’s  mantle. 
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Figure  1  Top  panel  shows  the  velocity  perturbations  with  respect  to  the  global  average 
at  150  km  depth  as  seen  in  S12WM13  for  the  Mediterranean  area.  The  bottom  panel  shows 
the  same  area,  but  for  the  model  S20U7L5. 
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Figure  2  Vertical  profiles  of  velocity  perturbation  with  respect  to  the  global  average  for 
two  locations  in  Figure  1  associated  with  large  velocity  perturbations,  West  Africa  and  the 
Red  Sea  Rift.  In  the  West  African  profile,  fast  velocities  continue  to  500  km  depth,  while 
the  Red  Sea  Rift  profile  shows  slow  velocities  down  .to  400  km  depth. 


S— velocity  (km/sec)  S— velocity  (km/sec) 


Figure  3  Vertical  velocity  profiles  for  two  locations  chosen  to  coincide  with  the  regions 
anah-zed  by  Grand  and  Helmberger  (1984)  in  the  development  of  the  models  SNA  and  TNA. 
Velocities  for  PREM,  SNA,  TNA,  and  S20U7L5  axe  shown  with  different  line  styles. 
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1.  Baltic  shield  (60N,  30E) 

2.  Northeast  Siberia  (SON,  120E) 

3.  Northwest  Africa  (20N,  10W) 

4.  Western  Australia  (30S,  120E) 

5.  Guyana  shield  (0,  SOW) 

S.  Canadian  shield  (SON,  90  W) 


S— velocity  (km/sec) 

Figure  4  Vertical  velocity  profiles  beneath  six  points  located  within  stable  continental 
interiors.  Also  shown,  with  the  thick  solid  line,  is  the  SNA  structure.  PREM  is  shown  with 
the  thin  solid  line. 
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